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LONG-TERM  GOALS 

The  long  term  goal  of  this  interdisciplinary  research  program  continues  to  be  the  demonstration  of 
Bayesian  Hierarchical  Model  (BHM)  utility  in  several  aspects  of  operational  ocean  forecasting. 

The  specific  goals  in  the  current  phase  of  the  research  are:  1)  publications  of  results  from  ensemble 
ocean  forecast  experiments  driven  by  the  surface  vector  wind  (SVW)  BHM;  b)  consolidating  impacts 
results  in  MFS  reforecast  experiments  for  the  time-dependent  error-covariance  BHM;  and  c)  running 
the  first  multi-model  and  multi-parameter  super-ensemble  BHM  experiments. 

OBJECTIVES 

Research  objectives  leading  to  the  publication  of  manuscripts  regarding  SVW-BHM  include: 

1)  constructing  3  appendices  for  Milliff  et  al  (2009)  to:  a)  demonstrate  a  systematic  approach  to 
(future)  process  model  development;  b)  document,  in  probability  model  notation,  the  complete 
SVW-BHM,  as  well  as  expressions  for  the  full  conditional  distributions;  and  c)  document  the 
SVW-BHM  hyperprior  specifications; 

2)  Re-writing  the  text  and  updating  figures  for  Bonazzi  et  al.  (2009);  and 
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3)  incorporating  co-author  final  edits  for  Milliff  et  al.  (2009)  and  Bonazzi  et  al.  (2009). 

Research  objectives  leading  to  consolidation  of  results  for  the  time-dependent  error-covariance  BHM 
include: 

1)  supplying  an  anomaly-only  data  stage  version  of  the  error-covariance  BHM  to  MFS  for 
reforecast  experiments;  and 

2)  interpreting  reforecast  results,  and  iterating  with  MFS  for  future  reforecast  experimental  design. 

Research  objectives  leading  to  the  first  multi-model  and  multi-parameter  super-enesmble  BHM  runs 
include: 

1)  evaluating  preliminary  experiments  based  on  a  Levantine  Intermediate  Water  (LIW)  formation 
rate  target  process; 

2)  re-organizing  the  target  process  to  focus  on  temperature  (7)  and  salinity  (S)  profile  evolution 
at  2  locations  in  the  region  of  LIW  formation  (Rhodes  Gyre);  and 

3)  providing  T(z,t)  and  S(z,t)  data  files  (.mat)  to  Mark  Berliner  for  preliminary  modelling. 

APPROACH 

Appendices  and  updates  for  the  BHM-SVW  papers: 

The  sample  BHM  model  equation  and  the  expression  for  the  full  conditional  distribution  for  the  zonal 
wind  component are  provided  here,  where  «  is  an  n-vector,  spanning  all  n  locations  on  the  target 
grid.  The  process  model  form  (derived  from  the  Rayleigh  Friction  Equations)  is: 


where  there  are  terms  in  the  Normal  distribution  specification  (right  hand  side)  for  pressure  gradient 
effects  and  a  scale-selective  nested  wavelet  summation  that  proved  useful  in  a  process  model  for 
tropical  winds  in  Wikle  et  al.  (2001).  In  keeping  with  a  key  advantage  of  BHM,  the  uncertainty  in  the 
zonal  wind  distribution  is  modeled  in  the  variance-covariance  term.  If  we  let  Hu  =  ajjDy  +  ai^2Dx  then 
the  full  conditional  distribution  is: 


(D^  -  -  hjiu) 


As  noted  above,  this  can  be  shown  to  be  a  Gaussian  Distribution  Ut\  ■  ~N(Ab,A),  where 


Expressions  of  these  forms  for  data  stage  and  process  model  terms,  and  the  implied  full  conditional 


distributions  of  BHM-SVW  are  included  in  Appendix  2  of  Milliff  et  al.  (2009).  The  hyperpriors  for 
terms  that  are  specified  (i.e.  rather  than  endowed  with  distributions)  in  the  full  conditional  distributions 
are  included  in  Appendix  3. 

The  companion  paper,  Bonazzi  et  al.  (2009)  was  updated  and  re-written  by  Prof  Nadia  Pinardi  during 
her  annual  summer  visit  to  Boulder,  CO  (CoRA  offices)  this  August.  Both  papers  are  undergoing  final 
edits  by  co-authors  before  submission  for  publication  this  month. 

Iterating  Time-Dependent  Error  Covariance  BHM  based  on  impacts: 

A  new  round  of  reforecast  experiments  for  the  Gulf  of  Lyon  region  (MFS  region  3)  have  been  run  to 
test  time-dependent  error  covariance  matrices  that  are  the  posterior  means  of  the  Error-Covariance 
BHM.  Earler  experiments  based  on  a  version  of  the  time-dependent  error  covariance  BHM  using  the 
quality  controlled  misfit  d  and  anomaly  q  data  provided  by  MFS  did  not  lead  to  significant 
improvements  in  comparisons  with  Argo  data.  The  comparisons  were  noisy,  suggesting  that  very 
localized  and  short  timescale  variability  was  influencing  the  time-dependent  error  covariance  in  ways 
that  were  not  useful  to  the  forecast  system. 

The  new  version  of  the  error  covariance  BHM  uses  only  q  information  in  the  data  stage.  The  q  inputs 
are  more  smoothly  varying  in  time  than  d,  and  they  exhibit  a  clear  seasonal  cycle.  Time  series  of  error 
covariance  matrices  were  provided  to  MFS  for  '^q-only"  reforecast  experiments. 

Identifying  T(z,t)>  S(z,t)  Target  Processes  for  the  Super-Ensemble  BHM 

The  LIW  formation  rate  time  series  from  analysis  (SYS3a2)  and  simulations  from  two  MFS  ocean 
models  (OPA  and  NEMO)  did  not  lead  to  informative  posterior  distributions  in  preliminary 
experiments  with  the  super-ensemble  BHM,  as  designed  and  run  by  Prof.  Mark  Berliner  (i.e.  as  in 
Berliner  and  Kim,  2008).  There  are  ambiguities  in  the  definition  of  LIW,  and  the  model  processes 
leading  to  formation,  that  are  larger  than  the  signal  and/or  model  biases.  As  such,  the  super-ensemble 
BHM  posteriors  were  dominated  by  the  model  uncertainties. 

During  Prof.  Nadia  Pinardi's  visit  to  CoRA  in  August,  she,  Berliner  and  Milliff  worked  with  Dr. 
Jeremiah  Brown  (NWRA/CoRA)  to  identify  two  locations  in  the  Rhodes  Gyre  region  where  time 
series  of  T  and  S  profiles  from  SYS3a2,  OPA  and  NEMO  could  be  used  as  a  new  target  process  for  the 
super-ensemble  BHM.  Time  series  of  T(z,t)  and  S(z,t)  have  been  provided  to  Prof.  Berliner  in 
preparation  for  new  super-ensemble  BHM  experiments.  The  thinking  is  that  inherent  uncertainties  in 
simulations  and  analysis  of  T(z,t)  and  S(z,t)  will  not  be  large  enough  to  dominate  posterior  distributions 
in  the  new  experiments  to  be  run  this  month. 


WORK  COMPLETED  AND  RESULTS 

Research  program  progress  and  plans  were  reviewed  for  ONR  program  managers  in  a  presentation  by 
Milliff  at  the  ONR  PI  Meeting  in  Chicago  in  June  2009. 

Publishing  BHM-SVW  development  and  forecast  impacts 

The  careful  production  of  the  BHM-SVW  equations  and  expressions  for  full  conditional  distributions, 
and  the  presentation  of  hyperprior  specifications,  are  essential  steps  in  preparing  BHM-SVW  for 
technology  transfer  to  operational  ocean  forecast  applications  at  MFS  and  elsewhere  (e.g.  NRL  SSC 
Ocean  Modeling  Group).  Together,  these  appendices  comprise  the  necessary  information  to  construct 
(and  adapt  as  needed)  the  BHM-SVW  used  in  our  demonstration  project  with  MFS. 

Bonazzi  et  al.  (2009)  demonstrates  impacts  on  MFS  forecasts  of  uncertainty  quantification  based  on 


realistic  spread  in  surface  wind  forcing  from  BHM-SVW.  Figure  1  provides  a  pictorial  sense  of  the 
enhanced  spread  in  ocean  model  response.  The  left  hand  panels  depicts  the  SSH  field  at  day  10  of 
ocean  forecasts  (i.e.  day  lOF)  driven  by  members  of  the  ECMWF  ensemble  forecast  winds,  and  the 
right  hand  panel  depicts  the  SSH  fields  on  lOF  as  forced  by  3  realizations  of  BHM-SVW.  Superposed 
on  each  SSH  field  are  6  trajectories  from  idealized  drifters  released  in  each  forecast  at  day  IF.  The 
mesoscale  response  in  the  BHM-SVW  case  exhibits  greater  spread,  as  do  the  simulated  trajectories. 
Ensemble  spreads  from  traditional  methodologies  notoriously  under-estimate  uncertainties  and  often 
do  not  exhibit  realistic  variability.  This  failing  is  apparently  obviated  by  BHM-SVW.  Quantification 
of  this  improvement  is  documented  in  Bonazzi  et  al  (2009). 


a)  ECMWF  1 


b)  BHM  1 


39.6N 

39.3N 


d)  BHM  2 

I _ I  - - 


c)  ECMWF  2 

I _ I  - 


e)  ECMWF  3 


f)  BHM  3 

LI— I  ~  - - 


15 

10 

5 

0 


-10 


Figure  1.  Sea-surface  height  anomaly  at  day  10  of  ensemble  ocean  forecasts  driven  by  (left)  members  of  the  ECMWF 
ensemble  forecast  winds  and  (right)  members  of  the  BHM-SVW  posterior  distribution.  Color  contours  indicate  SSH  in  cm, 
and  black  curves  mark  the  trajectories  of  6  idealized  drifters  released  in  the  forecast  fields  at  the  forecast  initial  time.  SSH 
variability  on  mesoscales,  and  spread  in  the  trajectories  are  larger  in  the  BHM-SVW  ensemble.  These  results  are  quantified 
in  Bonazzi  et  al.  (2009).  They  are  consistent  with  the  localization  and  mesoscale  focus  in  ocean  forecast  initial  conditions 
documented  in  Milliff  et  al.  (2009). 


Impacts  of  the  q-only  time-dependent  error  covariance  BHM 

Figure  2  depicts  a  comparison  of  root-mean  square  differences  in  sea-level  anomaly  (SLA)  with 
respect  to  an  operational  merged  sea-level  product  based  on  altimeter  data.  Three  forecast  experiments 
are  compared  in  MFS  region  3,  covering  the  same  period,  differing  only  in  the  form  of  the  error 
covariance  matrix  used  in  the  data  assimilation  step.  The  green  line  represents  the  operational  system 
wherein  the  error  covariance  matrix  is  fixed  in  time.  The  black  line  is  for  a  reforecast  using  time- 
dependent  error  covariance  matrices  from  the  BHM  based  on  both  d  and  q  data  stages,  and  the  yellow 
line  is  for  the  recent  q-only  experiments.  Surprisingly,  time-dependent  error  covariance  does  not  seem 
to  improve  the  comparisons  with  SLA  in  either  case.  Further  analyses  (not  shown)  demonstrate  that 
much  of  the  uncertainty  in  SLA  is  coming  from  the  error  covariances  for  salinity  (i.e.  as  opposed  to 
temperature).  New  designs  for  the  error  covariance  BHM  are  being  designed  with  these  findings  in 
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mind. 

Figure  2.  Sea-level  anomaly  (SLA)  root-mean-square  differences  with  respect  to  an  operational  satellite  merged  product  for 
MFS  model  grid  points  in  region  3  (Gulf  of  Lyon)  for  the  period  March-May  2005.  Three  reforecast  experiments  were  run 
over  the  period,  varying  only  the  form  of  the  error  covariance  matrix  used  in  the  data  assimilation.  The  green  line 
represents  the  operational  system  wherein  the  error  covariance  matrix  is  fixed  in  time.  The  black  line  represents  a  time- 
dependent  error  covariance  that  is  the  posterior  mean  output  of  the  error  covariance  BHM  based  on  misfit  (rf)  data  and 
anomaly  {q)  data.  The  yellow  line  corresponds  to  time-dependent  error  covariances  from  the  BHM  based  on  only  the  q 
data.  The  q-only  reforecast  experiment  only  ran  to  1  May  2005. 

Target  field  for  the  super-ensemble  BHM 

Figure  3  shows  sample  T(z,t)  summaries  for  a  point  on  the  southwestern  flank  of  the  Rhodes  Gyre  in 
the  Eastern  Mediterranean  Sea.  The  summary  period  shown  runs  from  1  February  through  31  March 
2007.  T(z,t)  and  S(z,t)  summaries  for  two  points  (the  other  location  is  on  the  northeastern  flank  of  the 
Rhodes  Gyre)  have  been  created  from  the  analysis  (1  realization  for  each  location)  and  from  1 1 
simulations  each  of  the  OPA  and  NEMO  development  models  at  MFS.  The  OPA  and  NEMO 
simulations  are  forced  by  the  ECMWF  analysis  and  10  realizations  of  the  BHM-SVW.  Similar 
summaries  are  provided  for  the  same  months  in  2005  and  2006.  The  super-ensemble  BHM  will  be 
used  to  hindcast  the  2008  T(z,t)  and  S(z,t)  summaries. 
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Figure  3.  T(z,t)  for  a  point  at  26.875°E,  33.5°N  (i.e.  on  the  southwestern  flank  of  the  Rhodes  Gyre)  for  the  period  1 
February  through  31  March  2007.  The  left  hand  panel  depicts  T(z,t)  from  the  MFS  analysis  (SYS3a2)  and  the  right  hand 
panel  is  from  an  OPA  model  simulation  driven  by  one  realization  from  BHM-SVW.  An  intermediate  water  formation  event 
is  evident  in  the  analysis  (left)  around  day  7.  It  is  absent  from  the  OPA  simulation  (right).  Evidence  of  a  mesoscale 
anticyclonic  eddy  appears  in  the  analysis  (left)  near  day  50.  The  mesoscale  variability  in  the  OPA  simulation  is  much 
reduced  vs.  the  analysis. 

IMPACTS  AND  APPLICATIONS 

Wikle  and  Milliff  were  invited  speakers  at  a  National  Research  Council  workshop  on  “Uncertainty 
Management  in  Remote  Sensing  of  Climate  Data”  in  December  2008.  The  BHM-SVW  was  the  topic 
of  Milliff  s  talk. 

Wikle,  Milliff  and  Dr.  Manuel  Fiadeiro  (ONR  Pgm  Mgr)  are  hosting  a  session  on  Probabilistic 
Modelling  in  Oceanography  at  the  AGU  2010  Ocean  Sciences  Meeting  in  Portland,  OR.  Berliner  and 
Pinardi  are  invited  speakers. 

RELATED  PROJECTS 

NSF-C433P-RFM-CORA- 1433: 

"Collaborative  Research:  Estimating  Ecosystem  Model  Uncertainties  in  Pan-Regional  Syntheses  and 
Climate  Change  Impacts  on  Coastal  Domains  of  the  North  Pacific  Ocean".  Co-PIs  include  Wikle  and 
Milliff;  Advisory  Council  includes  Berliner.  Period  of  Performance  September  2008-  August  2011. 

ONR  PO  (in  process): 

“Bayesian  Hierarchical  Model  Characterization  of  Model  Error  in  Ocean  Data  Assimilation  and 
Forecasts”  Co-PIs  include  Berliner,  Wikle  and  Milliff.  Period  of  Performance  to  be  determined 
(4yrs). 
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